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Abstract 

We present a method for numerically building a vortex knot state in the superfluid wave-function 
of a Bose-Einstein condensate. We integrate in time the governing Gross-Pitaevskii equation to 
determine evolution and stability of the two (topologically) simplest vortex knots which can be 
wrapped over a torus. We find that the velocity of a vortex knot depends on the ratio of poloidal 
and toroidal radius: for smaller ratio, the knot travels faster. Finally, we show how unstable vortex 
knots break up into vortex rings. 
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47.32.C- Vortex dynamics (fluid flow); 
03.75.Lm Vortices in Bose-Einstein condensation; 
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I. INTRODUCTION 



In 1867, following the works of Helmholtz on vortices and of Riemann on abelian functions, 
Lord Kelvin modelled atoms as knotted vortex tubes in the ether [Ij, effectively giving birth 
to knot theory [2j. This discipline has fascinated mathematicians and physicists since. More 
recently, knots have been the studied in different branches of physics, ranging from excitable 
media [3j, and classical field theory [4J, to optics [5^ 6j and liquid-crystal colloids [7J. 

Knots in superfluids are identified with closed vortex lines, regions of fluid around which 
the circulation assumes non-zero (quantized) value. Vortex rings have been studied experi- 
mentally in superfluid liquid helium |8lf9] and in Bose-Einstein condensates [lOj. Numerical 
simulations have revealed that superfluid turbulence contains linked vortex lines [llj, but, 
to the best of our knowledge, individual vortices with non-trivial topology have never been 
observed directly. To shed light onto this problem, energy, motion and stability of vortex 
knots have been examined theoretically and numerically using the classical theory of thin- 
cored vortex filaments. In this approach, the governing incompressible Euler dynamics is 
expressed by the Biot-Savart law or by its local induction approximation (LIA) [T2[ [T3] . 
Under certain conditions, it is found that some vortex knots are structurally stable, that is 
to say they travel without breaking up for distance larger than their own diameters. 

In superfluid helium, the validity of the classical theory of thin-core vortex filaments is 
based on the large separation of scales between the vortex core radius Qq (approximately 
10~^ cm in ^He and 10~^ cm in ^He-B) and the typical distance £ between vortices. In 
turbulence experiments, £ ^ 10~^ to 10~^ cm; the last value is also the typical diameter 
of experimental vortex rings [9j. The situation is very different in atomic Bose-Einstein 
condensates, where £ is only few times larger than a^. In this context, the Gross-Pit aevskii 
Equation (GPE) is clearly a more realistic model [14], particularly at very low temperatures, 
as thermal effects can be neglected. 

The advantage of the GPE is that it does not need the cut-off parameter required by the 
classical vortex filament theory to de-singularise the Biot-Savart integral The second 
advantage is that the GPE naturally describes vortex reconnections [16], which must be 
implemented algorithmically in the Biot-Savart model. Any prediction about the stability 
and the break-up of a vortex structure which is not orders of magnitude bigger than ao is 
therefore more reliable if obtained using the GPE. The third advantage of searching for vortex 
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knots in a Bose-Einstein condensate is that direct images of individual vortex structures are 
possible without the use of tracer particles which will certainly disturb these structures. 
The disadvantage is that atomic condensates are small, thus the motion of these structures 
will be affected by the boundaries and by the non-uniformity density of the background 
condensate. Before investigating these effects, however, it is essential to establish whether 
vortex knot solutions of the governing GPE exist, and, if they do, if they are sufficiently 
stable. This is the limited aim which we set in this work. We stress that we do not aim 
to propose a mechanism to experimentally create vortex knots in condensates, but only to 
study the possible existence and stability of these solutions of the GPE. We shall see that 
even setting up numerically a topologically non-trivial structure in the wave-function is not 
a minor task; indeed, to the best of our knowledge, this is the first time it has been done. 

The manuscript is organized as the following. Section [TT] explains how to create an 
elementary vortex knot in the initial conditions of the condensate wave-function. Section 



TM deals with the analysis of the dynamical properties of vortex knots. Section IV describes 



the break-up of vortex knots. Finally, the conclusions are in Section V 



II. VORTEX KNOT INITIAL CONDITIONS 

We consider the GPE written in the following dimensionless form 

2idtij-\/^^+\ij\^ij = 0, (1) 

where no external confining potential is present. The characteristic length scale of pertur- 
bations of the uniform condensate, called healing lengthy is defined as 

e = ^, where (p) = ^ [ (2) 
\/{p) ^ Jv 



>> 

is the mean density of the condensate. Besides the energy, the GPE conserves the total 
number of particles, therefore is a conserved quantity too. Without loss of generality, we 
choose to deal with a system that has an unperturbed density (the density field at infinity) 
equal to unity, and assume that perturbations are localized in a small region of the sample. 
In this hypothesis, ^ 2::^ 1 in our units. 

We now explain how to numerically build a vortex knot. First we construct a vortex. 
Consider the two dimensional plane sOz. A stable vortex is a hole (zero value) in the density 
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FIG. 1. (Colours online). The density field p around an axisymmetric two-dimensional vortex. 
Radial distances R are in units of the healing length ^. 

field, around which the phase of the wave-function changes by ±27r. A sufficiently accurate 
description of a two-dimensional vortex centered in the origin of the sOz plane is given by 
the wave- function ^2d{s-, z) = ^/p(R) e^^^^'^\ where R = V^5^~^F^, 



and the coefficients ai = 11/32, a2 = 11/384, hi = 1/3, and 62 = 11/384 arise from a second 
order Fade approximation [17j. Fig. [l] shows how the density field behaves around the 
axisymmetric vortex centre. It is clear from the plot that the vortex core is of the order of 
the healing length, and the bulk value of the density p = 1 is recovered at larger distances. 

We now come back to vortex knots in a three-dimensional system. We define a knot as 
a closed curve over a torus, characterized by the toroidal radius Rq and the poloidal radius 
Ri. More precisely, a closed curve Tn,m on the torus is determined by counting the number 
of toroidal wraps, n, and the number of poloidal wraps, m. For example, the curves 7i,i and 
72,2 describe respectively the unknot (the simple vortex ring) and two unlinked rings. The 
first topologically non-trivial curve is the trefoil, 72,3- In this work we shall focus on the two 
simplest knots, the trefoil 72,3 and its dual 73,2- 

A. The 7^,3 knot (trefoil) 



p{R) 



R^ (ai + a2R^) 
1 + + 62i?4 



(3) 



9{s^ z) = arctan 




The vortex line of a 72,3 knot lays on the torus as shown in Fig. [2j Any plane sOz passing 
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FIG. 2. (Colours online). Construction of the trefoil knot 72,3- 
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FIG. 3. (Colours online). Positions of the four vortices on = used to construct the wave- 
function of the trefoil knot 72,3- 

through the z axis intercepts the curve 72,3 at four different points, which correspond to four 
two-dimensional point vortices on the plane sOz. The positions of these two-dimensional 
vortices vary with respect to the choice of the plane sOz] in other words, these positions 
are functions of the angle variable introduced in Fig. [2] For example, the vortex positions 
of the wave-function for the angle = are shown in Fig. [3| By construction, these point 
vortices are located on the circumference defined by the intersection of the plane with the 
torus, and rotate on it following a particular function /((/>). To assure continuity of the vortex 
line and to describe the trefoil knot, the function /((/)) must have the form /((/)) = 30/2, 
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with (/) G [0,7r). 

We are now ready to write the three-dimensional wave-function which describes the tre- 
foil knot 72,3- In the approximation that the healing length ^ is much smaller than the 
inter-vortex distance, the two-dimensional wave-function in the plane sOz is given by the 
superposition (multiplication) of the wave-function ^^2D of each two-dimensional vortex cen- 
tred in the correct position, where the opposite circulation is obtained by applying the 
complex conjugation operator (...)*. Thus, the three-dimensional wave- function results in 



'02,3(^, z) = ^^2D {s{x, y) - Rq- Ri cos [a{x, y)],z- Ri sin y)]} 

X ^^2D {'^(x, y) — Rq — Ri cos [Qf(x, y) + tt] ^ z — Ri sin [a{x^ y) + ^]} 

X ^^2D {^{^^ y) + Ro — Ri cos y) + 7r/2] , z + i?i sin y) + 7r/2]} 

X ^^2^ {^{^^ y) + Rq — Ri cos [a{x, y) + 37r/2] ,z + Ri sin + 37r/2]} , 

(4) 

with ^) = sgn(x)^x2 + 7/^, where sgn(...) is the sign function, and 7/) = 3/2 arctan 



B. The 75,2 knot 

The technique used to define the wave-function of the trefoil knot can be extended to any 
other knot built on a torus. The 75,2 knot can be represented on the torus as shown in Fig. 
[4} In this case the generic plane sOz intersects the knot in six points, where the centers are 
function g{(f)) of the angle (f) and rotate around the circumference defined by the plane and 
the torus intersection. An example of the configuration for = is shown in Fig. [5} The 
function g{e) is g{e) = 29/3, with (p G [0,7r). 

Again, using the two-dimensional vortex description ^^2^, in the limit of inter- vortex 
distance much greater than the healing length the three-dimensional wave-function of a 
72,3 knot is 
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FIG. 4. (Colours online). Construction of the trefoil knot 73,2- 
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FIG. 5. (Colours online). Positions of the six vortices on = used to construct the wave- function 
of the trefoil knot 73,2- 
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V^3,2(^, z) = ^^2D {s{x, y) - Rq- Ri COS [a{x, y)],z- Ri sin [a{x, y)]} 

X ^^2D {s{x^ y) — Rq — Ri COS ^) + 27r/3] , z — i?i sin [q{(x, ^) + 27r/3]} 
X ^^2D {'^(x, y) — Rq — Ri cos + 47r/3] , z — i?i sin + 47r/3]} 

X ^^2D {'^(^5 y) + Rq + Ri COS [Qf(x, ^)] , z + Ri sin [q{(x, ^)]} 
X ^^2D {'^(^^ ?/) + ^0 + ^1 COS [a{x, y) + 27r/3] , z + i?i sin [a{x, y) + 27r/3]} 
X ^^2^ {^{^^ y) + Rq + Ri cos [a{x, y) + 47r/3] ,z + Ri sin [q{(x, + 47r/3]} , 

(5) 

with s{x^ y) = sgn(x)^x^ + y'^ and a{x^ y) = 2/3arctan (y/x). 

III. VORTEX KNOT DYNAMICS 

To study the dynamic and stabihty of the knots 72,3 and 75,2 with different geometries we 
need to find a compromise between the accessible numerical resolution and the box size: we 
need to resolve small scales near the vortex cores and, at the same time, minimize the finite 
size (boundary) effects. We recall that the parameters which identify our vortex knots, the 
toroidal and poloidal radii Rq and are expressed in units of the healing length ^. 

We chose to uniformly discretize physical space using a Cartesian grid with steps Ax = 
= Az = 0.5^ spanning over the knot ratios Ri/Rq = 1/10,1/5,2/5,3/5. We expect 
vortex knots to behave similarly to vortex rings, that is to say we expect that they travel 
along the direction of the torus axis of symmetry (the z-axis). Taking our computational 
constraints into account, we use 192 x 192 x 512 grid points {Lx = Ly = 96 ^ and = 256 (^) 
and the toroidal radius Rq = 20^. This choice allows us to have a minimum value of 
Ri = 2^ (when Ri/Rq = 1/10), acceptable to observe the small inter- vortex interactions, 
and a maximum knot size of 2{Rq + Ri) = 64(^ (when Ri/Rq = 3/5) which gives tolerable 
boundary effects. Table [T] summarizes the simulation parameters. 

In order to let the knot travel for the maximum distance in the z direction, at the start 
of the calculation {t = 0) the vortex knot is centered at the point [Lx/2^ Ly/2^ 2{R0 + RI)]. 
With this choice, the knot can propagate for a distance of 3 to 53/11 4.8 times its 
maximum size, before hitting the opposite side of the computational domain corresponding 
to z = Lz. The CPE is integrated in time using a split-step method with anti-periodic 
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(reflective) boundary conditions. The integration time step is At = 0.02 smaller that the 
fastest linear period — 0.032. This value allows us to conserves the initial energy and mass 
up to 3% and 1% respectively in all the simulations. Details on the numerical algorithm can 
be found in Ref. [18]. 

Fig. [6] and [7] show the iso-surfaces of the density field corresponding to the threshold 
value pth = 0.2 at the initial conditions and at successive times for the 72,3 and 73,2 knots 
respectively (unstable knots will be discussed in the next section). As expected, vortex 
knots move along the z-direction (the axis of symmetry of the torus), but also twist around 
it. Qualitatively, vortex knots with small knot ratio Ri/Rq are fast and stable, as they 
propagate along the z-direction without breaking. During the evolution, Kelvin waves [8j 
appear; such waves are visible at the last stages of cases (a), (b), (e), and (f). 

In order to quantify the evolution of vortex knots and compare one knot with others, we 
define the knot center of mass rcM = {xcm-> Vcm^ ^cm) as 

_ /,,rg(p,>-|Vf) dV 



where H{...) is the Heaviside step function. Fig. [8] and [9] show the z component zcm of 
the knot center of mass (shifted with respect to the initial position) for the 72,3 and 75,2 
cases respectively. In both cases, knots with smaller knot ratio Ri/Rq move faster, are more 
stable and propagate for longer distances before breaking up (a filled symbol at the end of 
each curve marks the break-up point). 

The z component of the velocity of a vortex knot is estimated by evaluating Vz{t) = 



[zcM{t + t) — zcM{t)] /t (where r = 4 for numerical convenience). Fig. 10 and Fig. O 



case 


knot ratio 


max size 


min size 


stable 


72,3 — 75,2 


Ri/Ro 


2(i?o + Ri) 


2Ri 


Yes/No 


(a) - (e) 


1/10 


44^ 


4C 


Y — Y 


(b) - (f) 


1/5 


48^ 


8^ 


Y — Y 


(c) - (g) 


2/5 


56^ 


16^ 


N — N 


(d) - (h) 


3/5 


64 C 


24 C 


N — N 



TABLE I. Vortex knot parameters of 72,3 and 73,2 used in the simulations. 
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(a) 



(b) 



(c) 



(d) 



FIG. 6. (Colours online). Iso-surfaces of the density field at the threshold level pth = 0.2 for 72,3 
knots of various knot ratios Ri/Rq (see Table [l]). Snapshots at times t = 0, 400, 800, 1200. Unstable 
knots are not shown. 



show Vz{t) measured in units of the vortex ring velocity [13 



^ring 



(R) 



n K 



In 



8R 

T 



0.615 



(7) 



having quantum number n = 1 and radius R = Rq (note that in our non-dimensional 
system the quantum of circulation is = 27r). It s apparent that vortex knots move with 
approximately constant z velocity before either breaking up or reaching the boundary of the 
computational domain, where the interaction with the image slows them down. 

It is instructive to analyze the mean and the standard deviation of the vortex knots' 
velocities measured in the constant-velocity regimes. The results, expressed in units of 

Three conclusions can be drawn from this figure: 



^rmg 



(i?o), are shown in Fig. 
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1. 72,3 knots are slower than 73,2 knots with the same knot ratio. This is physically 
expected as the velocity field of torus knots at large distance is similar to the velocity 
field of vortex rings with multiple circulation: 72,3 corresponds to circulation of 2tv 
and 73,2 to 3k. According to Equation ([7|), the velocity is directly proportional to the 
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FIG. 7. (Colours online). Iso-surfaces of the density field at the threshold level pth — 0.2 for 75,2 
knots of various knot ratios Ri/Rq (see Table [l]). Snapshots at times to t = 0,400,800. Unstable 
knots are not shown. 
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FIG. 8. (Colours online). Position along the z axis (in units of the healing length) of the centre of 
mass of 75,3 knots of various knot ratios Ri/Rq as a function of time. The filled circles denote the 
position where vortex knots break up. The horizontal lines denote respectively the distance where 
boundary effects become non-negligible and the finite system size along z. 
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Finite system size 
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FIG. 9. (Colours online). As in Fig.lslbut for 75,2- 
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FIG. 10. (Colours online). Velocity component Vz of 72,3 versus time of knots with various knot 
ratios Ri/Rq before destroying (filled points). Velocities are expressed in units of vortex ring 
velocity Q with quantum number n = 1 and radius R = Rq. 

circulation, and so 72,3 knots should be slower than 75,2 ones. However, this simple 
consideration does not apply well to knots because we would have expected, for the 
small knot ratio tested {Ri/Rq = 1/10) a scaled velocity of '^z/'^ring(^o) — ^ 
^^/^ring(^o) — 2 ^^.3 ^^.2 respectively, and this is not the case. 



2. The z- velocity component scales with the knot ratio and can be parametrised as 

Kq I Kq 



(8) 



where ^n,m and 5^,rn are coefficients which refer to the generic torus knot Tn,m- Values 



of Aji^m and Bn^m for the knots 72,3 and 73,2 are reported in Fig. 12, It is interesting 
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FIG. 11. (Colours online). Velocity component Vz of 73,2 knots having different knot ratio Ri/Rq 
before destroying (filled points). Velocities are expressed in units of ring velocity Q having quan- 
tum number n = 1 and radius R = Rq. 




R1/R0 



FIG. 12. (Colours online). Averaged velocity components Vz of 72,3 and 73,2 vortex knots with 
various knot ratios Ri/Rq. Velocities are expressed in units of vortex ring velocity Q with quantum 
number n = 1 and radius R = Rq; error-bars correspond to one standard deviation. 



to observe that yl2,3 — ^3,2- 



3. The z-velocity component of unstable knots (i.e. knots that decay before reaching the 
computational boundaries) is less or similar to Vj;[j^g{Ro)' On the contrary, knots that 
remain stable within the computational domain are characterized by > '^^•^^^(^o)- 
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FIG. 13. (Colours online). Three successive snapshots showing how the 72,3 vortex knot with knot 
ratio Ri/Rq = 2/5 breaks up into two vortex rings. Here we plot two perspectives (up and side to 
the vortex propagation) of the iso-surfaces of the density field corresponding to the threshold level 
Pth = 0.2. 

IV. THE BREAKING OF A KNOT 

In our simulations we have observed that some vortex knots break up into topologically 
simpler objects. 

We first analyze the unstable 72,3 knots: these are knots corresponding to knot ratios 
Ri/Ro = 2/5,3/5. In Fig. 13 we show three snapshots of the decay of the 72,3 knot 
with ratio Ri/Rq = 2/5. It is apparent that the knot breaks into two vortex rings via 
three simultaneous self-reconnection events (see in particular the snapshot corresponding to 
t = 224). The decay of the vortex knot 72,3 with ratio Ri/Rq = 3/5, not shown here, is 
similar. 



On the contrary, 75,2 vortex knots break in a different manner. As shown in Fig. 14 



the vortex knot 75,2 with knot ratio Ri/Rq = 3/5 first decays in one vortex ring and 
two linked vortex rings via two simultaneous self-reconnection events (snapshot at time 
t = 192). Subsequently, the small free ring escapes from the other rings, which undergo 
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t=180 



t=192 



t=204 



t=416 



t=428 



t=440 



FIG. 14. (Colours online). Two sequences of three successive snapshots showing how the 73,2 vortex 
knot with knot ratio Ri/Rq = 3/5 breaks up into three vortex rings. Here we plot two perspectives 
(up and side to the vortex propagation) of the iso-surfaces of the density field corresponding to the 
threshold level pth = 0.2. 



two simultaneous reconnection events that create two unlinked vortex rings (snapshot at 
time t = 428). The last step is remarkable: there is no apparent reason why two linked 
vortex knots should in principle unlink into two vortex rings (by making two simultaneous 
reconnection events) without forming a single ring (by one reconnection event). 



The 75,2 vortex knot with ratio Ri/Rq = 2/5 qualitatively decays in the same way, 
producing a set of three unlinked vortex rings, but the steps are quite different. In the first 
step, a free vortex ring and two linked vortex rings are again produced. However, the free 
ring, which is initially located behind the two linked vortex rings, is smaller and faster than 
the other rings. As a consequence, it reconnects with the two linked vortex rings, as shown 
in Fig. [15] (snapshot t = 752). At this point the reformed knot breaks up, undergoing 
the same sequence previously described for the 75,2 with ratio Rq/Ri = 3/5 case, and the 
outcome is a set of three vortex rings (snapshot at time t = 952). Note that, in the last 
snapshot, the first knot (in the sense of the position) has split into two smaller vortex rings 
via a self-reconnection event which is probably consequence of its Kelvin waves oscillations. 
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t=552 t=752 t=952 

FIG. 15. (Colours online). Three successive snapshots showing how the 73,2 knot with knot ratio 
Ri/Rq = 2/5 breaks into three vortex rings. Here we plot two perspectives (up and side to the 
vortex propagation) of the iso-surfaces of the density field corresponding to the threshold pth = 0.2. 



V. CONCLUSIONS 



We have numerically analyzed the existence and stability of vortex knots in the GPE 
model of a condensate. We have proposed a novel numerical technique for creating ab initio 
vortex knots in the wave-function of the condensate. In particular, we have focussed our 
numerical computations on the two simplest (in the topological sense) vortex knots, 72,3 and 
the 73,2- We have analyzed the evolution and the stability of such knots with respect to the 
knot ratio Ri/Rq. We have found that a knot can be unstable, i.e. it breaks up into simple 
rings during the propagation, or stable within our computational domain. Our numerical 
experiments clearly show that a small knot ratio {Ri/Rq = 1/10, 1/5) increases the stability, 
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whereas a large knot ratio {Ri/Rq = 2/5,3/5) decreases it, in agremeent with [T2]. 

We have found that vortex knots propagate essentiaUy as vortex rings. We have measured 
the vortex knot velocities along the torus symmetry axis and shown that the velocity depends 
linearly on the knot ratio for both 72,3 and 73,2- 

Finally, we have studied the details of the break up of vortex knots. Although we do 
not have a theoretical explanation for the break up, we have observed evidences of generic 
breaking behavior: 72,3 knots always break into two vortex rings via a three simultaneous 
self-reconnection event, whereas. 72,3 knots first decay into three vortex rings via two simul- 
taneous self-reconnect ions which create a free ring and two linked rings, then undergo two 
simultaneous reconnections which split the resulting link. 

We believe that our work opens up new interesting problems in the field of fiuid topology 
applied to superfiuids and Bose-Einstein condensates. The natural developments of our 
study will be a theoretical investigation of the stability of vortex knots, and an experimental 
study of the creation of a knotted initial condition in an atomic condensate. 
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